#include "header.h"

static int 	myid, np;
static char 	line[80];
fftw_complex	a;
double		Ntot, ptime;
double		d1, d2, d3, c1, c2, c3, c4;

int main(int argc, char **argv)
{
  ctrl_str	ctrl;
  work_str	work;
  int 		p, flag;
  MPI_Status    status;
  MPI_Request   request;

  d1 = 1 / ( 2 - pow(2, 1/3.) );
  d3 = d1;
  d2 = 1 - (d1+d3);

  c1 = 0.5 * d1;
  c2 = 0.5 - c1;
  c3 = 0.5 - c1;
  c4 = c1;
  
  MPI_Init(&argc, &argv);
  //fftw_mpi_init();
  MPI_Comm_size(COMM, &np);
  MPI_Comm_rank(COMM, &myid);

  load_params(&argc, &argv, &ctrl);
  arr3d_init(&ctrl, &work);
  init_io(&ctrl);
  init_ffts(&ctrl, &work);
  init_timers(&ctrl);
  init_splitstep();
  set_ic(&ctrl, &work);
  if (ctrl.gFlag == 0) ft_potential();
  
  if (ctrl.rs_flag == 1) load_data(&work, ctrl.rname);
  //if (myid == 0) printf("time multiplier %f + (%f)I\n",creal(step(&ctrl)),cimag(step(&ctrl)));
  if (ctrl.m_order == 2) lin_step(&work, 0.5*step(&ctrl)); // for 2-nd order
  if (ctrl.m_order == 4) nnl_step(&work, c1*step(&ctrl));   // for 4-th order
  while ( stop_evolve() ){
    time_out(&work, &ctrl);
    if (ctrl.m_order == 2) s_order(&ctrl, &work);
    if (ctrl.m_order == 4) f_order(&ctrl, &work);
    adv_time();
    save_temporal(&work,"phase_t.out"); 
    Ntot = particles(&work);
    //if (myid == 0) printf("\tnum particles %.15e\n", Ntot);
    //time_out(&work, &ctrl); 
  } 
  save_data(&work, "final.bin");
  fftw_mpi_cleanup();  
  MPI_Finalize();
  return 0;

}
